function impedance_from_a(CONSTS, plot_data)

    phi = 0;
    m_max = 20; %20
    k0 = CONSTS.k0;
    eps = CONSTS.eps;
    eta = CONSTS.eta;
    eps_a = CONSTS.eps_a;
    c = CONSTS.c;
    Z0 = 4*pi/c;
    coef_sgs_to_si   = 9e+11;
    
    a_low_order  = 0;
    a_high_order = 5;
    num_points   = 1000;
    a_vec = [0.5; 1; 2; 3; 4; 5; 6; 7; 8; 9; 10]; %linspace(a_low_order,a_high_order,num_points)'; %

    Z_vec = zeros(size(a_vec));
    Z_an_vec = zeros(size(a_vec));
    Z_odn_pl_vec = zeros(size(a_vec));
    for i = 1:size(a_vec,1)
        a = a_vec(i);
        d = a/50;
%         fprintf('a=%d\n', a);

        LOCAL_CONSTS = CONSTS;
        LOCAL_CONSTS.a = a;
        LOCAL_CONSTS.d = d;
        LOCAL_CONSTS = update_CONSTS(LOCAL_CONSTS);
        fprintf('a=%d\n', a);
        fprintf('d=%d\n', d);
        
        % input impedance obtained numerically
        I_sum_0 = func_I_sum(phi, m_max, LOCAL_CONSTS);
        Z_vec(i) = coef_sgs_to_si ./ I_sum_0;

        % input impedance obtained analytically
%         eps_ef = ((2i*eps_a)/((1+(eps_a)^2/(abs(eps*eta)))*(sqrt(abs(eps*eta))))) - 4;
        eps_ef = (eps_a^2+abs(eps*eta))/((2*sqrt(abs(eps*eta))-1i*eps_a));
%         h1 = k0*(abs(eps*eta))^(1/4)*((1+1i)/sqrt(2*eps_ef));
        h1 = k0*sqrt(eps_ef)*(1-1i)/sqrt(2);
        I_sum_analytical = -((1i*pi*h1)/(Z0*k0*log(4*a/d)))*((cos((pi-abs(phi)).*h1.*a))./(sin(pi*h1*a)));
        Z_an_vec(i) = coef_sgs_to_si/I_sum_analytical;
        
        % input impedance in the case of a homogeneous magnetoplasma
        h = k0*(abs(eps*eta))^(1/4)*((1-1i)/sqrt(2)); 
        I_sum_odn = -((1i*pi*h)/(Z0*k0*log(4*a/d)))*((cos((pi-abs(phi)).*h.*a))./(sin(pi*h*a)));
        Z_odn_pl_vec(i) = coef_sgs_to_si/I_sum_odn;      
        
        if(false)
            fprintf('a     =  %e  cm \n', a);
            fprintf('Z     =  %e  + i %e  Om \n', real(Z_vec(i)), ...
                                                  imag(Z_vec(i)));
            fprintf('Z_odn =  i %e  Om \n',       imag(Z_odn_pl_vec(i)));
        end
    end

    if(plot_data)
        figure; plot(a_vec, real(Z_vec), 'b.-', ...
                     a_vec, imag(Z_vec), 'bx-', ...
                     a_vec, real(Z_an_vec), 'r-', ...
                     a_vec, imag(Z_an_vec), 'r-', ...
                     a_vec, real(Z_odn_pl_vec), 'm-', ...
                     a_vec, imag(Z_odn_pl_vec), 'm-'); 
        title('Z(a)'); 
        xlabel('a, cm'); ylabel('R, X');
        legend('Re Z', 'Im Z', 'Re Z_{an}', 'Im Z_{an}', 'Re Z_0', 'Im Z_0', ...
            'location', 'SouthEast');
    end

end